%% Figures 1-6 & 8

clc;
clear;

%% Loading data

load('Data\Special_Transactions_Pseudo.mat');             
Specials.VolRate=Specials.Rate.*Specials.Volume; 
load('Data\GC_Transactions_Pseudo.mat');
GC.VolRate=GC.Rate.*GC.Volume;
StoreFolder='Figures\EpsFiles';

% download data on the DFR (see e.g.: https://www.ecb.europa.eu/stats/policy_and_exchange_rates/key_ecb_interest_rates/html/index.en.html) and replace code below:
DFR = table();%this creates pseudo-data
t1 = datetime(2010,1,1);%this creates pseudo-data
t2 = datetime(2019,1,1);%this creates pseudo-data
DFR.refDate = [t1:t2]';%this creates pseudo-data
DFR.DFR = 0.5.*ones(1,length(DFR.refDate))';%this creates pseudo-data

clear t1 t2

%% Figure 1: Different market turnovers

% Step 1: Data for Figure 1

% Repo volumes
opts = spreadsheetImportOptions("NumVariables", 3);
% opts.Sheet = "Sheet1";
opts.DataRange = "A2:C21";
opts.VariableNames = ["quarter", "GC", "Total"];
opts.VariableTypes = ["datetime", "double", "double"];
opts = setvaropts(opts, "quarter", "InputFormat", "");
RepoQuarterlyMarketVolumes = readtable("Data\RepoQuarterlyMarketVolumes_Pseudo.xlsx", opts, "UseExcel", false);

clear opts

% Download data on unsecured trading volumes 
% From the Euro Money Market Survey (EMMS) until 2015, and replace code below:
UnsecCumQ2015 = table();%this creates pseudo-data
t1 = datetime(2010,1,1);%this creates pseudo-data
t2 = datetime(2015,1,1);%this creates pseudo-data
UnsecCumQ2015.Date = [t1:calmonths(12):t2]';%this creates pseudo-data
UnsecCumQ2015.QuarterlyVol = 10000000.*ones(1,length(UnsecCumQ2015.Date))';%this creates pseudo-data

clear t1 t2

% From the Money Market Statistical Reporting (MMSR) after 2017, when available, and replace code below:
UnsecCumQA2017 = table();%this creates pseudo-data
t1 = datetime(2017,1,1);%this creates pseudo-data
t2 = datetime(2019,1,1);%this creates pseudo-data
UnsecCumQA2017.Date = [t1:calmonths(3):t2]';%this creates pseudo-data
UnsecCumQA2017.QuarterlyVol = 6000000.*ones(1,length(UnsecCumQA2017.Date))';%this creates pseudo-data

clear t1 t2

% Download data total cumulative QE purchases under the Public Sector Purchase Programme (PSPP) and volumes of the ECB’s main refinancing
% operations (MRO) and targeted longer-term refinancing operations (TLTRO), and replace code below:
VolumesCB = table();%this creates pseudo-data
t1 = datetime(2010,1,1);%this creates pseudo-data
t2 = datetime(2019,1,1);%this creates pseudo-data
VolumesCB.Date = [t1:calmonths(3):t2]';%this creates pseudo-data
VolumesCB.TotalPSPP = 150000.*ones(1,length(VolumesCB.Date))';%this creates pseudo-data
VolumesCB.MROTLTRO = 100000.*ones(1,length(VolumesCB.Date))';%this creates pseudo-data

clear t1 t2

% Step 2: Create figure

Fig = figure;
plot(VolumesCB.Date,VolumesCB.TotalPSPP,'Color',[0.71 0.7 0.7],'LineWidth',2)
hold on;
plot(VolumesCB.Date,VolumesCB.MROTLTRO,'Color',[0.93 0.79 0],'LineWidth',2)
hold on;
plot(RepoQuarterlyMarketVolumes.quarter,RepoQuarterlyMarketVolumes.Total,'g','LineWidth',2)
hold on;
plot(RepoQuarterlyMarketVolumes.quarter,RepoQuarterlyMarketVolumes.GC,'r','LineWidth',2)
hold on;
plot(UnsecCumQ2015.Date,UnsecCumQ2015.QuarterlyVol,':b','LineWidth',2)
hold on;
plot(UnsecCumQA2017.Date,UnsecCumQA2017.QuarterlyVol,'b','LineWidth',2)
xlim([datetime('01.01.2010','InputFormat','dd.MM.yyyy') datetime('01.01.2019','InputFormat','dd.MM.yyyy')])
legend('PSPP','MRO + TLTRO','Secured Market - GC&Special','Secured Market - GC','Unsecured EMMS','Unsecured MMSR','Location','NorthWest')

clear Fig RepoQuarterlyMarketVolumes UnsecCumQ2015 UnsecCumQA2017 VolumesCB
 
%% Figure 2: Repo rate development

% Step 1: Data for Figure 2

% Average GC rate
DayGC = varfun(@nansum, GC(GC.Country=='DE',:), 'InputVariables', {'VolRate' 'Volume'}, 'GroupingVariables', {'ECBAccessLender' 'refDate'});
DayGC.VolWeightAvgRate = DayGC.nansum_VolRate./DayGC.nansum_Volume;
DayGC = removevars(DayGC, {'GroupCount' 'nansum_VolRate','nansum_Volume'});
DayGC = sortrows(DayGC,'refDate','ascend');

% Average special rates for QE eligible and noneligible assets during PSPP
SpecialsFig2 = Specials;
SpecialsFig2.Eligible = 0.*ones(height(SpecialsFig2),1);
SpecialsFig2.Eligible(SpecialsFig2.refDate > datetime('08.03.2015','InputFormat','dd.MM.yyyy')) = SpecialsFig2.HypEligible(SpecialsFig2.refDate > datetime('08.03.2015','InputFormat','dd.MM.yyyy'));
DaySpecials = varfun(@nansum, SpecialsFig2(SpecialsFig2.Country=='DE',:), 'InputVariables', {'VolRate' 'Volume'}, 'GroupingVariables', {'Eligible' 'refDate'});
DaySpecials.VolWeightAvgRate = DaySpecials.nansum_VolRate./DaySpecials.nansum_Volume;
DaySpecials = removevars(DaySpecials, {'GroupCount' 'nansum_VolRate','nansum_Volume'});
DaySpecials = unstack(DaySpecials,'VolWeightAvgRate','Eligible','NewDataVariableNames',{'NonEligible','Eligible'});

% Step 2: Create figure

Fig = figure;
plot(DaySpecials.refDate,DaySpecials.NonEligible)
hold on;
plot(DaySpecials.refDate,DaySpecials.Eligible);
hold on;
plot(DayGC.refDate,DayGC.VolWeightAvgRate);
xlim([datetime('01.01.2010','InputFormat','dd.MM.yyyy') datetime('01.01.2019','InputFormat','dd.MM.yyyy')])
hold on;
plot(DFR.refDate, DFR.DFR);
yline(0);
xline([datetime('08.03.2015','InputFormat','dd.MM.yyyy')],'-',{'Begin QE'},'LabelHorizontalAlignment','Center')
legend('show');
ylim([-2 2])
xlim([datetime('01.01.2010','InputFormat','dd.MM.yyyy') datetime('01.01.2019','InputFormat','dd.MM.yyyy')])
legend('Noneligibile','Eligible','GC','Deposit Rate')

clear DayGC DaySpecials Fig SpecialsFig2
    
%% Figure 3: Repo market and deposit facility volumes
% Figure 3a: General collateral trading volumes

% Step 1: Data for Figure 3a

[G,refDate,ECBAccessLender] = findgroups(GC.refDate,GC.ECBAccessLender);
Vol = splitapply(@nansum,GC.Volume,G);
TradingVol = unstack(table(refDate,ECBAccessLender,Vol),'Vol','ECBAccessLender','NewDataVariableNames',{'NoAccessVol','AccessVol'});
TradingVol.TotalVol = TradingVol.NoAccessVol + TradingVol.AccessVol;
TradingVol.NoAccessShare = TradingVol.NoAccessVol./TradingVol.TotalVol;

% Step 2: Create figure

Fig = figure;
area(TradingVol.refDate,TradingVol.TotalVol);
hold on;
area(TradingVol.refDate,TradingVol.NoAccessVol);
xlim([datetime('01.01.2010','InputFormat','dd.MM.yyyy') datetime('01.01.2019','InputFormat','dd.MM.yyyy')]);
legend('Lender Access ','Lender Nonaccess');

clear ECBAccessLender Fig G refDate TradingVol Vol

% Figure 3b: Deposit facility volume and spread between GC rate and DFR

% Step 1: Data for Figure 3b

DayGC = varfun(@nansum, GC, 'InputVariables', {'VolRate' 'Volume'}, 'GroupingVariables', {'refDate'});
DayGC.VolWeightAvgRate = DayGC.nansum_VolRate./DayGC.nansum_Volume;
DayGC = removevars(DayGC, {'GroupCount' 'nansum_VolRate','nansum_Volume'});
DayGC = sortrows(DayGC,'refDate','ascend');    
DayGC = outerjoin(DayGC,DFR,'Type','left','MergeKeys',1,'LeftKey','refDate','RightKey','refDate');   
DayGC.Diff = DayGC.VolWeightAvgRate - DayGC.DFR;

% download data on the volumes placed at the DFR (e.g. from Bloomberg) and replace code below:
DepFacilityUsage = table();%this creates pseudo-data
t1 = datetime(2010,1,1);%this creates pseudo-data
t2 = datetime(2019,1,1);%this creates pseudo-data
DepFacilityUsage.refDate = [t1:caldays(7):t2]';%this creates pseudo-data
PseudoData= 150000000000.*[0.5:(1/length(DepFacilityUsage.refDate)):1.5]';%this creates pseudo-data
DepFacilityUsage.DepFacilityVol = PseudoData(1:end-1); %use pseudo-data
% Step 2: Create figure

Fig = figure;
yyaxis right;
h(2)=plot(DepFacilityUsage.refDate,DepFacilityUsage.DepFacilityVol);
ax = gca;
ax.YColor = [0, 0, 0];
yyaxis left;
h(1)=plot(DayGC.refDate(1:end),DayGC.Diff);
ylim([-0.5 0.5])
h(3)=yline(0);
legend(h([1 2]),'GC-Deposit Rate (lhs)','Deposit Facility Volume (rhs)')
ax = gca;
ax.YColor = [0, 0, 0];
xlim([datetime('01.01.2010','InputFormat','dd.MM.yyyy') datetime('01.01.2019','InputFormat','dd.MM.yyyy')])

clear ax DayGC DepFacilityUsage Fig h t1 t2

%% Figure 4: Market order share and excess rate frequency

% Step 1: Data for Figure 4

GCFig4 = outerjoin(GC,DFR,'Type','left','MergeKeys',1,'LeftKey','refDate','RightKey','refDate');
    DayGC = varfun(@nansum, GCFig4, 'InputVariables', {'VolRate' 'Volume'}, 'GroupingVariables', {'refDate' 'Term' 'Basket'});
    DayGC.VolWeightAvgRate = DayGC.nansum_VolRate./DayGC.nansum_Volume;
    DayGC = removevars(DayGC, {'GroupCount' 'nansum_VolRate','nansum_Volume'});
GCFig4 = outerjoin(GCFig4,DayGC,'Type','left','Keys',{'refDate' 'Basket' 'Term'},'MergeKeys',1);
GCFig4.IsPositiveER = (GCFig4.Rate - GCFig4.VolWeightAvgRate) > 0;
GCFig4.IsMarketOrderByLender=GCFig4.Aggressor=="Lender";

DayGCFig4 = varfun(@nanmean, GCFig4(GCFig4.Term=="ON",:), 'InputVariables', {'IsPositiveER' 'IsMarketOrderByLender'}, 'GroupingVariables', {'ECBAccessLender'});
DayGCFig4.Name = categorical(DayGCFig4.ECBAccessLender);
DayGCFig4.Name(DayGCFig4.ECBAccessLender==1) = "Access";
DayGCFig4.Name(DayGCFig4.ECBAccessLender==0) = "Nonaccess";
DayGCFig4.Name = removecats(DayGCFig4.Name);

% Step 2: Create figure

Fig = figure;
grayColor=[.7 .7 .7]; 
subplot(1,2,1)
bar(DayGCFig4.Name,DayGCFig4.nanmean_IsMarketOrderByLender,'FaceColor',grayColor)
ylim([0 1])
title('Market Order Share')
subplot(1,2,2)
bar(DayGCFig4.Name,DayGCFig4.nanmean_IsPositiveER,'FaceColor',grayColor)
ylim([0 1])
title('Lending Rate')

clear DayGC DayGCFig4 Fig GCFig4 grayColor

%% Figure 5: Trading share of access banks

% Step 1: Data for Figure 5

GCFig5 = GC;
GCFig5.refMonth = dateshift(GCFig5.refDate,'start','month');

GCFig5.German = GCFig5.Basket=='GC_Germany';%this uses pseudo-basket name
    MonthGCDE = varfun(@nansum, GCFig5(GCFig5.German==1,:), 'InputVariables', {'Volume'}, 'GroupingVariables', {'ECBAccessLender' 'refMonth'});
    MonthGCDE = removevars(MonthGCDE, 'GroupCount');
GraphDE = unstack(MonthGCDE,{'nansum_Volume'},{'ECBAccessLender'},'NewDataVariableNames',{'NoAccessVol','AccessVol'});
GraphDE.TotalVol = nansum([GraphDE.NoAccessVol,GraphDE.AccessVol],2);
GraphDE.AccessShare = GraphDE.AccessVol./GraphDE.TotalVol;
GraphDE = sortrows(GraphDE,'refMonth','ascend');

GCFig5.Core = GCFig5.Country=='DE'|GCFig5.Country=='FR'|GCFig5.Country=='BE'|GCFig5.Country=='NL'|GCFig5.Country=='FI'|GCFig5.Country=='AT';
    MonthGCCore = varfun(@nansum, GCFig5(GCFig5.Core==1 & GCFig5.German==0,:), 'InputVariables', {'Volume'}, 'GroupingVariables', {'ECBAccessLender' 'refMonth'});
    MonthGCCore = removevars(MonthGCCore, 'GroupCount');
GraphCore = unstack(MonthGCCore,{'nansum_Volume'},{'ECBAccessLender'},'NewDataVariableNames',{'NoAccessVol','AccessVol'});
GraphCore.TotalVol = nansum([GraphCore.NoAccessVol,GraphCore.AccessVol],2);
GraphCore.AccessShare = GraphCore.AccessVol./GraphCore.TotalVol;
GraphCore = sortrows(GraphCore,'refMonth','ascend');

    MonthGCAll = varfun(@nansum, GCFig5(GCFig5.Core==0 & GCFig5.German==0,:), 'InputVariables', {'Volume'}, 'GroupingVariables', {'ECBAccessLender' 'refMonth'});
    MonthGCAll = removevars(MonthGCAll, 'GroupCount');
GraphAll = unstack(MonthGCAll,{'nansum_Volume'},{'ECBAccessLender'},'NewDataVariableNames',{'NoAccessVol','AccessVol'});
GraphAll.TotalVol = nansum([GraphAll.NoAccessVol,GraphAll.AccessVol],2);
GraphAll.AccessShare = GraphAll.AccessVol./GraphAll.TotalVol;
GraphAll = sortrows(GraphAll,'refMonth','ascend');

% Step 2: Create figure

Fig=figure;
plot(GraphDE.refMonth,GraphDE.AccessShare)
hold on
plot(GraphCore.refMonth,GraphCore.AccessShare)
hold on
plot(GraphAll.refMonth,GraphAll.AccessShare)
hold on
legend('Collateral-driven','Funding-driven, more expensive','Funding-driven, less expensive','location','southwest')

clear Fig GCFig5 GraphAll GraphCore GraphDE MonthGCAll MonthGCCore MonthGCDE

%% Figure 6: Special collateral repo market
% Figure 6a: Special collateral trading volume

% Step 1: Data for Figure 6a

DaySpecials = varfun(@nansum, Specials, 'InputVariables', {'Volume' 'VolRate'}, 'GroupingVariables', {'refDate' 'HypEligible' });
DaySpecials.VolWeightAvgRate = DaySpecials.nansum_VolRate./DaySpecials.nansum_Volume;
DaySpecials = removevars(DaySpecials, {'GroupCount' 'nansum_VolRate'});
DaySpecials_Volume = unstack(DaySpecials(:,{'refDate' 'HypEligible' 'nansum_Volume'}),'nansum_Volume','HypEligible','NewDataVariableNames',{'Noneligible','Eligible'});
DaySpecials_Volume.Total = DaySpecials_Volume.Noneligible + DaySpecials_Volume.Eligible;

% Step 2: Create figure

Fig = figure;
area(DaySpecials_Volume.refDate,DaySpecials_Volume.Total);
hold on;
area(DaySpecials_Volume.refDate,DaySpecials_Volume.Eligible);
xlim([datetime('01.01.2010','InputFormat','dd.MM.yyyy') datetime('01.01.2019','InputFormat','dd.MM.yyyy')])
legend('Noneligible','Eligible');

clear DaySpecials_Volume Fig 

% Figure 6b: Spread between (hypothetically) eligible and noneligible assets

% Step 1: Data for Figure 6b

DaySpecials_Rate = unstack(DaySpecials(:,{'refDate' 'HypEligible' 'VolWeightAvgRate'}),'VolWeightAvgRate','HypEligible','NewDataVariableNames',{'Noneligible','Eligible'});
DaySpecials_Rate.RateDiff = DaySpecials_Rate.Noneligible - DaySpecials_Rate.Eligible; 

% Step 2: Create figure

Fig = figure;
plot(DaySpecials_Rate.refDate,DaySpecials_Rate.RateDiff);
yline(0);
xline([datetime('08.03.2015','InputFormat','dd.MM.yyyy')],'-',{'Begin QE'},'LabelHorizontalAlignment','Center');
xline([datetime('10.11.2015','InputFormat','dd.MM.yyyy')],'-',{'Expansion'},'LabelHorizontalAlignment','left');
xline([datetime('13.01.2017','InputFormat','dd.MM.yyyy')],'-',{'Extension'},'LabelVerticalAlignment','bottom','LabelHorizontalAlignment','left');
xlim([datetime('01.01.2010','InputFormat','dd.MM.yyyy') datetime('01.01.2019','InputFormat','dd.MM.yyyy')]);

clear DaySpecials DaySpecials_Rate Fig

%% Figure 8: Repo market dispersion

% Step 1: Data for Figure 8 
% DK index general collateral segment

DayGC = varfun(@nansum, GC, 'InputVariables', {'Volume' 'VolRate'}, 'GroupingVariables', {'refDate'});
DayGC.VolWeightAvgRate = DayGC.nansum_VolRate./DayGC.nansum_Volume;
DayGC = removevars(DayGC, {'GroupCount' 'nansum_VolRate','nansum_Volume'});

DayGC_Access = varfun(@nansum, GC, 'InputVariables', {'Volume' 'VolRate'}, 'GroupingVariables', {'refDate' 'ECBAccessLender'});
DayGC_Access.VolWeightAvgRate = DayGC_Access.nansum_VolRate./DayGC_Access.nansum_Volume;
DayGC_Access = removevars(DayGC_Access, {'GroupCount' 'nansum_VolRate'});
DayGC_Access = unstack(DayGC_Access,{'VolWeightAvgRate' 'nansum_Volume'},{'ECBAccessLender'},'NewDataVariableNames',{'NoAccessVol','AccessVol'});

GraphGC = outerjoin(DayGC_Access,DayGC,'Type','left','MergeKeys',1,'LeftKey','refDate','RightKey','refDate');  
GraphGC.Dispersion = (GraphGC.nansum_Volume_AccessVol.*abs(GraphGC.VolWeightAvgRate_AccessVol - GraphGC.VolWeightAvgRate) + GraphGC.nansum_Volume_NoAccessVol.*abs(GraphGC.VolWeightAvgRate_NoAccessVol - GraphGC.VolWeightAvgRate))./(GraphGC.nansum_Volume_AccessVol + GraphGC.nansum_Volume_NoAccessVol);
GraphGC = outerjoin(GraphGC,DFR,'Type','left','MergeKeys',1,'LeftKey','refDate','RightKey','refDate');
GraphGC.IsBelowDFR = GraphGC.VolWeightAvgRate < GraphGC.DFR;

clear DayGC DayGC_Access

% DK index special segment

DaySpecials = varfun(@nansum, Specials, 'InputVariables', {'Volume' 'VolRate'}, 'GroupingVariables', {'refDate'});
DaySpecials.VolWeightAvgRate = DaySpecials.nansum_VolRate./DaySpecials.nansum_Volume;
DaySpecials = removevars(DaySpecials, {'GroupCount' 'nansum_VolRate','nansum_Volume'});

DaySpecials_Eligible = varfun(@nansum, Specials, 'InputVariables', {'Volume' 'VolRate'}, 'GroupingVariables', {'refDate' 'HypEligible'});
DaySpecials_Eligible.VolWeightAvgRate = DaySpecials_Eligible.nansum_VolRate./DaySpecials_Eligible.nansum_Volume;
DaySpecials_Eligible = removevars(DaySpecials_Eligible, {'GroupCount' 'nansum_VolRate'});
DaySpecials_Eligible = unstack(DaySpecials_Eligible,{'VolWeightAvgRate' 'nansum_Volume'},{'HypEligible'},'NewDataVariableNames',{'NonEligible','Eligible'});

GraphSpecial = outerjoin(DaySpecials_Eligible,DaySpecials,'Type','left','MergeKeys',1,'LeftKey','refDate','RightKey','refDate');  
GraphSpecial.Dispersion = (GraphSpecial.nansum_Volume_Eligible.*abs(GraphSpecial.VolWeightAvgRate_Eligible - GraphSpecial.VolWeightAvgRate) + GraphSpecial.nansum_Volume_NonEligible.*abs(GraphSpecial.VolWeightAvgRate_NonEligible - GraphSpecial.VolWeightAvgRate))./(GraphSpecial.nansum_Volume_Eligible + GraphSpecial.nansum_Volume_NonEligible);

clear DaySpecials DaySpecials_Eligible

% Step 2: Create figure
Fig = figure;

% Upper panel
subplot(2,1,1);
% To create gray-shaded areas
a = area(GraphGC.refDate,GraphGC.IsBelowDFR.*max(GraphGC.Dispersion));
a.FaceAlpha = 0.05;
a.EdgeAlpha = 0;
a.FaceColor = 'k';
yline(0);
hold on;
plot(GraphGC.refDate,GraphGC.Dispersion);
ylim([0 max(GraphGC.Dispersion)]);
xlim([datetime('01.01.2010','InputFormat','dd.MM.yyyy') datetime('01.01.2019','InputFormat','dd.MM.yyyy')]);
legend('GC<DFR','DK Index','Location','northwest');
title('DK index general collateral segment');

%lower panel
subplot(2,1,2);
plot(GraphSpecial.refDate,GraphSpecial.Dispersion);
xline([datetime('08.03.2015','InputFormat','dd.MM.yyyy')],'-',{'Begin QE'},'LabelHorizontalAlignment','left','Color','k');
xline([datetime('10.11.2015','InputFormat','dd.MM.yyyy')],'-',{'Expansion'},'LabelHorizontalAlignment','left','Color','k');
xline([datetime('13.01.2017','InputFormat','dd.MM.yyyy')],'-',{'Extension'},'LabelVerticalAlignment','bottom','LabelHorizontalAlignment','left','Color','k');
ylim([0 max(GraphSpecial.Dispersion)]);
xlim([datetime('01.01.2010','InputFormat','dd.MM.yyyy') datetime('01.01.2019','InputFormat','dd.MM.yyyy')]);
title('DK index special segment');

clear a Fig GraphGC GraphSpecial